| event | clusterid | feature | mfi | |
|---|---|---|---|---|
| Min. : 776 | Length:554130 | Length:554130 | Min. :-1.9451 | |
| 1st Qu.: 9171 | Class :character | Class :character | 1st Qu.: 0.4495 | |
| Median :17176 | Mode :character | Mode :character | Median : 2.2822 | |
| Mean :19285 | NA | NA | Mean : 2.3179 | |
| 3rd Qu.:28082 | NA | NA | 3rd Qu.: 4.1336 | |
| Max. :47312 | NA | NA | Max. : 6.8623 |
Background
In this analysis workflow we use spectral flow cytometry to simultaneously assess histone H3 post-translational modifications (H3-PTMs) at the single-cell level, in the context of cell cycle status and cell-type-specific markers. This workflow was implemented in BioProtocol Exchange by Golden et al., 2024 and contains supplementary figures and data analysis associated with the pre-print.
Information about biological material
iPSC-derived Neuronal Progenitor Cells (NPCs) were generated using the female control line SCTi003-A (Stem Cell Technologies). NPC induction was performed with the STEMdiff SMADi Neural Induction Kit, following the manufacturer’s instructions. All experiments were conducted with three biological replicates, defined in this context as NPCs induced in three independent wells. For further details on experimental design and data acquisition, please refer to Golden et al., 2024.
Information about initial data processing performed in OMIQ
| Data processing | Strategy | Additional Information |
|---|---|---|
| Scaling | Hyperbolic arcsine (arcsinh) | Finak et al., 2010 |
| Gating | Cell QC and manual gating for comparison | Golden et al., 2024 |
| Subsetting | Fully stained samples | Removal of single color controls and blanks to prevent noise in downstream analysis |
| Subsampling | Random downsampling. Equal number of cells available in the Live Cells filter | Braanker et al., 2021 Note: Density-dependent downsampling strategy should be tested to assess if there is loss of small cell populations due to random cell selection Qiu, 2017 |
| Dimension Reduction | PCA | all features |
| Dimension Reduction | opt-SNE | selected features |
| Clustering | FlowSOM | selected features |
Data analysis and visualization
The goal of this analysis is to assess H3-PTM signatures in NPCs in the context of their cell cycle status and lineage markers (PAX6 in this case).
I named the fluorescence intensity column as mfi and it is misleading. Note that this is not the median value but the fluorescence value of each feature in each cell.
| event | clusterid | feature | mfi | mfi_z.V1 | |
|---|---|---|---|---|---|
| Min. : 776 | Length:554130 | Length:554130 | Min. :-1.9451 | Min. :-2.2360456 | |
| 1st Qu.: 9171 | Class :character | Class :character | 1st Qu.: 0.4495 | 1st Qu.:-0.9800129 | |
| Median :17176 | Mode :character | Mode :character | Median : 2.2822 | Median :-0.0186857 | |
| Mean :19285 | NA | NA | Mean : 2.3179 | Mean : 0.0000000 | |
| 3rd Qu.:28082 | NA | NA | 3rd Qu.: 4.1336 | 3rd Qu.: 0.9523977 | |
| Max. :47312 | NA | NA | Max. : 6.8623 | Max. : 2.3836750 |
Z-score column was added to data
Goal: Assess H3-PTM values in the context of cell cycle status and markers such us cleaved-Caspase3 and PAX6. This section also aims to compare FlowSOM unsupervised clustering and manual gates.
In the context of this project (and most spectral data) the targets of the antibodies are referred to as features. Gates, whether manual or created through FlowSOM, are referred to as filters.
Initial visualization of all events within each cluster
This plot allows visualization of single events in each cluster. Not included in Methods paper
Adapted from Heatmap for time series
Show code
p <-ggplot(df,aes(event,feature,fill=mfi_z))+
geom_tile(color= "white",linewidth=0.000000001) +
scale_fill_viridis(name="Scaled Fluorescense",
alpha = 1,
begin = 0,
end = 1,
direction = 1,
discrete = FALSE,
option = "C",
aesthetics = "fill")
p <-p + facet_grid(~clusterid, scales = "free", space = "free",)
p <-p + scale_y_discrete(limits = rev(levels(df$feature))) # Reverse the order of the factor levels for y-axis
p <-p + scale_x_discrete(labels = NULL)
p <-p + theme_minimal(base_size = 8)
p <-p + labs(title= paste("Fluorescence intentisity values of each feature by fsom cluster | NPC"), x="Event", y="Feature")
p <-p + theme(legend.position = "bottom")+
theme(plot.title=element_text(size = 14))+
theme(axis.text.y=element_text(size=6)) +
theme(strip.background = element_rect(colour="white"))+
theme(plot.title=element_text(hjust=0))+
theme(axis.ticks=element_blank())+
theme(axis.text=element_text(size=7))+
theme(legend.title=element_text(size=8))+
theme(legend.text=element_text(size=6))+
theme(strip.text.x = element_text(angle = 90, hjust = 1)) + # Rotate the cluster id labels
removeGrid()#ggExtra
# you will want to expand your plot screen before this bit!
pHeatmap visualization of concatenated files
For this plot, we used heatmap data downloaded from OMIQ, which includes three concatenated biological replicates. To compare the resolution between FSOM clusters and manual gating, data will be plotted for both sets of filters.
| H3K14ac | phH3_S10 | H3K4me3 | Total.H3 | Caspase3 | PAX6 | H3K27ac | H3K9ac | FxCycle | |
|---|---|---|---|---|---|---|---|---|---|
| fsom_02 | 0.45 | 0.24 | 4.80 | 0.80 | 0.03 | 1.25 | 3.00 | 4.12 | 3.57 |
| fsom_09 | 0.78 | 0.54 | 5.39 | 1.68 | -0.38 | 3.25 | 4.67 | 5.35 | 4.20 |
| fsom_04 | 0.47 | 0.23 | 4.85 | 0.93 | 1.42 | 2.14 | 3.30 | 4.34 | 3.63 |
| fsom_01 | 0.53 | 0.24 | 5.10 | 1.23 | -0.15 | 1.24 | 3.74 | 4.52 | 3.78 |
| fsom_05 | 0.79 | 0.37 | 5.44 | 1.69 | -0.35 | 1.78 | 4.44 | 5.24 | 4.23 |
| fsom_08 | 0.55 | 0.23 | 5.09 | 1.18 | -0.15 | 2.81 | 3.97 | 4.82 | 3.99 |
| H3K14ac | phH3_S10 | H3K4me3 | Total.H3 | Caspase3 | PAX6 | H3K27ac | H3K9ac | FxCycle | |
|---|---|---|---|---|---|---|---|---|---|
| Apoptotic Cells | 0.16 | 0.21 | 4.05 | 0.30 | 4.70 | 0.90 | 2.12 | 2.13 | 3.52 |
| fsom_02 | 0.45 | 0.24 | 4.80 | 0.80 | 0.03 | 1.25 | 3.00 | 4.12 | 3.57 |
| Other/G0/G1 | 0.46 | 0.23 | 4.90 | 0.93 | -0.03 | 1.25 | 3.25 | 4.26 | 3.64 |
| fsom_09 | 0.78 | 0.54 | 5.39 | 1.68 | -0.38 | 3.25 | 4.67 | 5.35 | 4.20 |
| Other/G2 | 0.71 | 0.27 | 5.35 | 1.24 | -0.31 | 1.35 | 4.01 | 5.05 | 4.36 |
| Other/Mitosis | 1.74 | 5.42 | 5.55 | 1.31 | -0.38 | 1.43 | 4.10 | 4.78 | 4.55 |
Features MFI by fsom clusters
Hierarchical Clustering (hclust)
Initial assessment through hierarchical clustering
Show code
[1] "#324DA0FF" "#008EA8FF" "#69BBABFF" "#C9DDC9FF" "#FEFDBEFF" "#F4D481FF"
[7] "#EB9D00FF" "#DE5400FF" "#A51122FF"
Show code
# Sort dendogram nodes
row_dend = dendsort(hclust(dist(df_m)))
col_dend = dendsort(hclust(dist(t(df_m))))
# Plot heatmap
h1 <- Heatmap(df_m, name = "MFI",
row_title_rot = 0,
row_names_gp = gpar(fontsize = 8), # Set row name font size
column_names_gp = gpar(fontsize = 8), # Set column name font size
col = col_1,
rect_gp = gpar(col = "white", lwd = 0.5),
cluster_rows = row_dend,
cluster_columns = col_dend,
row_split = 2, column_split = 2,
row_gap = unit(0.5, "mm"), column_gap = unit(0.5, "mm"), border = TRUE,
column_title = "MFI Values Across NPC Clusters Identified by FlowSOM | hclust",
column_title_gp = gpar(fontsize = 8, fontface = "bold"),
cell_fun = function(j, i, x, y, width, height, fill) {
grid.text(sprintf("%.1f", df_m[i, j]), x, y, gp = gpar(fontsize = 5)) # Add text inside cells
})
h1quartz_off_screen
2
quartz_off_screen
2
Pearson Correlation
Useful to cluster data based on the correlation between features rather than their absolute differences. It’s particularly good for identifying clusters with similar expression patterns across different markers, regardless of absolute values.
Show code
[1] "#324DA0FF" "#008EA8FF" "#69BBABFF" "#C9DDC9FF" "#FEFDBEFF" "#F4D481FF"
[7] "#EB9D00FF" "#DE5400FF" "#A51122FF"
Show code
# Plot heatmap
# Create the heatmap with Pearson correlation as the distance metric
set.seed(135)
h2 <- Heatmap(df_m, name = "MFI",
clustering_distance_rows = "pearson", # Use Pearson correlation for row clustering
clustering_method_rows = "complete", # Specify the clustering method (e.g., complete linkage)
row_km = 4, row_km_repeats = 100,
row_title_rot = 0,
row_names_gp = gpar(fontsize = 9), # Set row name font size
column_names_gp = gpar(fontsize = 9), # Set column name font size
col = col_1,
row_gap = unit(1, "mm"), border = TRUE,
column_title = "MFI Values Across NPC Clusters Identified by FlowSOM | pearson",
column_title_gp = gpar(fontsize = 8, fontface = "bold"),
cell_fun = function(j, i, x, y, width, height, fill) {
grid.text(sprintf("%.1f", df_m[i, j]), x, y, gp = gpar(fontsize = 7)) # Add text inside cells
})
h2quartz_off_screen
2
quartz_off_screen
2
Features MFI in manual gates
Show code
[1] "#3C5941FF" "#3C5941FF" "#758A68FF" "#B4BD93FF" "#FBF2C4FF" "#E3C28DFF"
[7] "#D48E53FF" "#C7522BFF" "#C7522BFF"
Show code
# Sort dendogram nodes
row_dend = dendsort(hclust(dist(dfg_m)))
col_dend = dendsort(hclust(dist(t(dfg_m))))
# Plot heatmap
h1g <- Heatmap(dfg_m, name = "MFI",
col = col_fall,
rect_gp = gpar(col = "white", lwd = 0.5),
cluster_rows = row_dend,
cluster_columns = col_dend,
row_split = 7, column_split = 2,
row_gap = unit(1, "mm"), column_gap = unit(1, "mm"), border = TRUE,
column_title = "NPC MFI by manual gating | hclust")
h1gquartz_off_screen
2
quartz_off_screen
2
This data viz strategy aims to highlight the contributions of each feature to different filters (aka clusters).
Prepare data
Show code
# Prepare data
cbp <- data_long
# Convert relevant columns to factors
cbp <- cbp %>%
mutate(
feature = factor(feature),
clusterid = factor(clusterid),
event = factor(event),
mfi_shifted = mfi + abs(min(mfi)) + 1 # Shift MFI (not scaled) values to avoid negatives
)
summary_cbp <- summary(cbp)
kable(summary_cbp, format = "html", caption = "Summary of data containing shifted values")| event | clusterid | feature | mfi | mfi_shifted | |
|---|---|---|---|---|---|
| 2595 : 27 | fsom_06:182997 | Caspase3: 61570 | Min. :-1.9451 | Min. :1.000 | |
| 2598 : 27 | fsom_08:116640 | FxCycle : 61570 | 1st Qu.: 0.4495 | 1st Qu.:3.395 | |
| 2600 : 27 | fsom_01: 88902 | H3K14ac : 61570 | Median : 2.2822 | Median :5.227 | |
| 2602 : 27 | fsom_02: 83502 | H3K27ac : 61570 | Mean : 2.3179 | Mean :5.263 | |
| 2604 : 27 | fsom_09: 39312 | H3K4me3 : 61570 | 3rd Qu.: 4.1336 | 3rd Qu.:7.079 | |
| 2605 : 27 | fsom_05: 19179 | H3K9ac : 61570 | Max. : 6.8623 | Max. :9.807 | |
| (Other):553968 | (Other): 23598 | (Other) :184710 | NA | NA |
Features contributions to each cluster
To assess circular barplots for all clusters, please go to the output folder.